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ABSTRACT 

We present a new technique to quantify cluster-to-cluster variations in the observed 
present-day stellar mass functions of a large sample of star clusters. Our method 
quantifies these differences as a function of both the stellar mass and the total cluster 
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mass, and offers the advantage that it is insensitive to the precise functional form of the 
mass function. We applied our technique to data taken from the Advanced Camera for 
Surveys Survey for Globular Clusters, from which we obtained completeness-corrected 
1 stellar mass functions in the mass range 0.25-0.75 M Q for a sample of 27 clusters. The 

results of our observational analysis were then compared to Monte Carlo simulations 
for globular cluster evolution spanning a range of initial mass functions, total numbers 
of stars, concentrations, and virial radii. 

We show that the present-day mass functions of the clusters in our sample can 
be reproduced by assuming an universal initial mass function for all clusters, and 
\ that the cluster-to-cluster differences are consistent with what is expected from two- 

body relaxation. A more complete exploration of the initial cluster conditions will be 
needed in future studies to better constrain the precise functional form of the initial 
mass function. This study is a first step toward using our technique to constrain the 
dynamical histories of a large sample of old Galactic star clusters and, by extension, 
star formation in the early Universe. 

?H \ Key words: globular clusters: general - stellar dynamics - stars: statistics - methods: 

statistical - stars: formation - stars: low-mass. 



1 INTRODUCTION 

It is now thought that most, if not all, o f the stars in 
our Galaxy were born in star cluste rs (e.g. ILada fc Ladal 
1 19951 . 120031 : IMcKee fe Ostrikerj l2007l b And yet, there re- 
main several key details of the star formation process 
that are still not understood. Part of the problem lies 
in the fact that populations of young stars are typi- 
cally hidden by a dense veil of optically-thick gas and 
dust. This prevents the escape of most of the light pro- 
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duced by infant stars, and often renders these regions dif- 
ficult to observe (e.g. [Grenier. Casandiian Sz Terrier! 120051 ; 
ILada. Alves fe Lombardil l2007l b Most of these clusters are 
sparsely p opulated an d are of relatively low mass (M < 10 4 
M Q ) (e.g. lLada|[l985T ). They are also very young since clus- 
ters of suc h low mass are unlikely to survive for more than 
1 Gyr fe.g. lPortegies Zwart. McMillan fe Gieles!l201ol '). 

At the other end of the cluster mass spectrum, most 
massive star clusters (M > 10 4 Mq) in our Galaxy tend to 
be at least a few Gyrs old, and in many cases are nearly 
as old as the Universe itself (e.g. lHarrislll996. 2010 update! ; 
iDe Angeli et~aH 120051 ) . These clusters have the advantage 
that they are no longer obscured by the primordial gas 
from which they formed, however they are a dynamically 
active environment. As a result, the conditions present at 
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the time of their formation have now been largely erased 
fe.g.lPortegies Zwart et al.ll200ll ; iHurlev et al]|2005l ; iMurravl 
I2009T ). This presents a considerable challenge for studying 
star formation in the regime of cluster masses and metal- 
licities that characterize Milky Way globular clusters. This 
is unfortunate since these old star clusters contain the fos- 
sil record of a very early episode of star formation in the 
Universe, and are the only means of studying it locally in 
massive star clusters. 

One of the primary observational tests for star forma- 
tion theories is the stellar initial mass function (IMF). Cur- 
rent observational evidence suggests that the IMF is very 
similar in different regions o f our Galaxy, including the disk 
and y oung star clusters fe.g. lElmegreenlll999l ; lKroupa et al.l 
l201lf ). Howev er, this is still being debated throughout the 
litera ture (e.g.[Scalo 1998l; |Parravano. McKee fc Hollenbacbl 
l201ll ). Different star formation theories tend to predict 
different IMFs. These vary with the properties of the 
gas clouds from which the stars are bor n, including den- 
sity, temperature and comp ositio n (e.g. lElmegreenl [200lj 
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Given the sensitive nature of the observations, a large 
sample of IMFs spanning the entire range of cluster proper- 
ties exhibited by star clusters in the Milky Way, including 
total mass and chemical composition, has yet to be compiled. 
This is a sorely needed step in order to advance our under- 
standing of star formation by providing direct comparisons 
between observations and theoretical predictions. This is 
especially true of massive, metal-poor star clusters since we 
are particularly lacking observations of I MFs in this regime 
of clu st er masses and metallicites (e.g. McKee fc Ostrikerl 
120071 : IPorteeies Zwart. McMillan fc Gieled l2010l \ IrrT 
port ant steps in this direction were rec ently taken 



by De Marchi. Paresce fc Portegies Zwartl l|201(t ) and 



iPaust et al.l (|2010l ). who studied the present-day mass func- 
tions of a large sample of Galactic clusters and considered 
the effects of the cluster dynamics in modifying them from 
their primordial forms. 

For the very first time, the Advanced Camera for Sur- 
veys (ACS) Survey for Globular Clusters has provided pho- 
tometry for a large sample of Milky Way globular clusters 
(GCs) that reaches down to unprecedented faint magni- 
tudes. This offers a large sample of current stellar mass func- 
tions spanning the stellar mass range « 0.2 — 0.8 Mq. All of 
the clusters are massive and very old, with total masses and 
ages ra nging from fts 10 4 - 10 6 Mq and ~ 10-12 Gyrs, respec- 
tively (|Harris! Il996. 2010 update!; |Pe Angeli etafl l2005h . 
This has allowed significant time for their stellar mass func- 
tions to have been modified from their primordial forms 
due to both stellar evolution and stellar dynamics. How- 
ever, most of the processes responsible for this evolution are 
now largely understood. Therefore, in principle, it is possi- 
ble to use current observations of old star clusters together 
with theoretical models for their evolution to extrapolate 
backwards in time and indirectly probe their IMFs. 

For most of the life of a massive star cluster, two- 
body relaxation is the dominant physical mechanism driving 
its evolution (e.g. lHenonlll960l; ISpitzerlll987l ; iHeggie fc Hut] 
120031 ; iGieles, Heggie fc Zhadl201ll '). The term describes the 
cumulative effects of long-range gravitational interactions 
that occur between pairs of stars, which act to alter their or- 



bits within the cluster. This results in a phenomenon known 
as mass segregation, which is the tendency for heavier stars 
to accumulate in the central cluster regions and low-mass 
stars to be dispersed to wider orbits. This mechanism also 
causes stars to escape from their host cluster, with the prob- 
ability of ejection increasing with decreasing cluster mass. 
Therefore, two-body relaxation acts to slowly modify the 
distribution of stellar masses within clusters, and can cause 
very dynamically evolved clusters to appear severely de- 
pleted of their low-mass stars. Evidence in favour of this pro- 
cess having actually occurred in r eal star clusters has been 
reported by several authors (e.g. Ivon Hippel fc Saraiediixil 



ll99Sl : lDe Marchi. Paresce fc Portegies ZwarlJl2010h 

A number of theoretical studies have been conducted to 
learn how the evolution of the stellar mass function (MF) 
in GCs is affected by two-body relaxation, stellar evolu- 
tion, disc shocking, and tida l effects from the Galaxy (see 
iBaumgardt fc Makinol l|2003h for a detailed review). In the 
absence of these effects, we expect the MF to continually 
rise toward lower st ellar masses. By performin g a series of 
N-body simulations. I Vesperini fc Heggie! (|l997f ) showed that 
tidal effects from the Galaxy, disc shocking, and a higher 
initial central concentration all act to increase the rate of 
stellar evaporation, and accelerate the depletion of preferen- 
tially lo w-mass stars. These results we re confirmed and built 
upon bv lBaumgardt fc Makinol (|2003l 1 who showed that the 
depletion of low-mass stars can be sufficiently dramatic to 
change the sign of the slope of the MF at the low-mass 
end. Interestingly, th ese results were not supported by th e 
observational study of lDe Marchi. Paresce fc Pulonel l|2007h . 
These authors analyzed the MFs in a sample of 20 Galactic 
GCs, and found that the slope of the MF decreases with 
increasing central concentration. They argued that this con- 
tradicts what is expected from theory since two-body relax- 
ation is responsible both for increasing the central density 
and flattening the MF at the low-mass end. In an effort 
to explain this, they suggested that many of the clusters 
in their sample could be post-core collapse, and therefore 
had much higher central densi ti es in the past. Alternatively, 
iMarks. Kroupa fc Baumgardtj (1 20081 ) argued that this can 
be explained by residu al gas-expulsion from initially mass 
segregated clusters (e.g. lTutukov|[l97^ ) . and cautioned that 
unresolved binaries could also be contributing. 

Several theoretical studies have also been conducted 
to study the dynamical histories of indi vidual globular 
clusters (e.g. IHeggie fc Gierszl |2008L I2009D . For example, 
IZonoozi et al.l (2011) recently performed the first ever di- 
rect N-body simulations of a Milky Way (MW) GC over its 
entire lifetime. This was done for the distant GC Palomar 
14, which has an unusually low-density and large radius. The 
emphasis of this paper is to use the ensemble information of 
many GCs to learn about the universality of the IMF in old 
massive star clusters. Individual cases, in particular Pal 14, 
are often chosen for their peculiar characteristics and may 
not be representative of the bulk of the GCs in the MW. 

In this paper, we present a new technique to quan- 
tify cluster-to-cluster variations in the observed stellar mass 
functions of a large sample of clusters spanning a diverse 
range of properties. Our method offers the advantage that 
it is insensitive to the precise functional form of the MF. 
We have applied it to a sample of 27 MFs taken from th e 
ACS Survey for Globular Clusters |Saraiedini et all [20071 '). 
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Can the present-day MFs be explained by an universal IMF 
and stellar evaporation induced by two-body relaxation? Or 
are cluster-specific IMFs needed to reproduce the observed 
MFs? To address these questions, we compared the results of 
our observational analysis to 268 Monte Carlo simulations 
for GC evolution. The models spanned a range of initial 
masses, virial radii, central concentrations and IMFs. There- 
fore, by evolving all of these models to the current ages of 
the GCs in our sample and comparing the resulting MFs to 
the present-day observed ones, we have quantified the dy- 
namical evolution of the MFs in our observed sample. This 
has allowed us to take the first steps toward constraining 
both the exact functional forms of the IMFs of MW GCs 
and the conditions present at the time of their formation. 

In Section we present our sample of observed stellar 
mass functions and describe both our technique for ana- 
lyzing the observations and the models for globular cluster 
evolution to which they are compared. The results of our 
analysis of the ACS observations are presented in Section [3] 
along with an example comparison between the observations 
and the models. This example demonstrates how our method 
can be used to compare a large number of observed MFs to 
analogous samples of simulated MFs. Finally, we discuss in 
Section U the implications of our results for the conditions 
present in our observed clusters at the time of their forma- 
tion and the role played by two-body relaxation in modifying 
the stellar MF to its present-day form. 



2 METHOD 

In this section, we describe how we acquired our sample of 
mass functions from the ACS data, as well as the Monte 
Carlo simulations for globular cluster evolution used for 
comparison to the observations. 

2.1 The Data 

The data used in t his study was taken from the sample of 35 
MW GCs used in iLeigh. Sills fc Kniggel (|201ll ), which was 
in turn taken from th e ACS Survey for Globular Clusters 
l|Saraiedini et al 1120071) Fl The ACS Survey provides unprece- 
dented deep photometry in the F606W (« V) and F814W 
(~ I) filters that is nearly complete down to ~ 0.2 Mq. 
In other words, the colour-magnitude diagrams (CMDs) ex- 
tend reliably from the horizontal branch all the way down 
to about 7 magnitudes below the main-sequence turn-off 
(MSTO). A list of the GCs used in this study is shown 
in Table Q] along with their core radii (r c ), half-mass radii 
(r^), central luminosity densities (po), and central con centra- 
tion parameters (c) . These were taken directly from lHarrisI 
l|l996. 2010 update!) , with the exception of the core and half- 
mass radii. The latter quantities are given in parsecs and 
were calculated using the distance modulii and e xtinction 
corrections provided in lHarrisI l|l996. 2010 update]) . 

Each cluster was centred in the ACS field, which ex- 
tends out to several core radii from the cluster centre in 
most cases. Coordinates for the cluster centres were taken 

1 The data can be found at 

http: ~www.astro.ufl.edu/^ata/public_hstgc/| 



Table 1. List of globular clusters and their structural parameters 



Cluster 


Alternate 


r c 




P0 


c 


I T 

1JJ 


LL) 


(in pc) 


(in pc) 


(in L Q pc -3 ) 




104 


47 Tuc 


0.47 


4.11 


4.88 


2.07 


1261 




1.66 


3.22 


2.99 


1.16 


2298 




0.97 


3.08 


2.90 


1.38 


4147 




0.51 


2.70 


3.63 


1.83 


4590 


M 68 


1.73 


4.51 


2.57 


1.41 


5024 


M 53 


1.82 


6.80 


3.07 


1.72 


5272 


M 3 


1.10 


6.84 


3.57 


1.89 


5286 




0.95 


2.48 


4.10 


1.41 


5904 


M 5 


0.96 


3.85 


3.88 


1.73 


5927 




0.94 


2.46 


4.09 


1.60 


5986 




1.43 


2.97 


3.41 


1.23 


6093 


M 80 


0.44 


1.78 


4.79 


1.68 


6121 


M 4 


0.59 


2.18 


3.64 


1.65 


6171 


M 107 


1.04 


3.21 


3.08 


1.53 


6205 


M 13 


1.29 


3.51 


3.55 


1.53 


6218 


M 12 


1.11 


2.49 


3.23 


1.34 


6254 


M 10 


0.98 


2.49 


3.54 


1.38 


6304 




0.36 


2.43 


4.49 


1.80 


6341 


M 92 


0.63 


2.45 


4.30 


1.68 


6535 




0.71 


1.68 


2.34 


1.33 


6584 




1.02 


2.86 


3.33 


1.47 


6637 


M 69 


0.84 


2.15 


3.84 


1.38 


6779 


M 56 


1.21 


3.02 


3.28 


1.38 


6838 


M 71 


0.74 


1.96 


2.83 


1.15 


6934 




1.00 


3.14 


3.44 


1.53 


6981 


M 72 


2.28 


4.60 


2.38 


1.21 


7089 


M 2 


1.08 


3.56 


4.00 


1.59 



from lGoldsburv et all (|2010h . These authors found their cen- 
tres by fitting a series of ellipses to the density distributions 
within the inner 2' of the cluster centre, and computing an 
average value. 

2.2 Measuring the Stellar Mass Function 

First, we used the available photometry to obtain estimates 
for the masses of the stars in our sample. To do this , we fit 
theoretical isochrones taken from lDotter et all (|2007l ) to the 
CMDs of every cluster. Each isochrone was generated using 
the metallicity and age of the cluster, and fit to its CMD 
using the c orresponding distan ce modulus and extinction 
provided in iDotter et~aU l|2010l ). The MSTO was then de- 
fined using our isochrone fits by selecting the bluest point 
along the main-sequence (MS). 

We considered five stellar mass bins along the MS. 
These ranged from 0.25 - 0.75 M© in increments of 0.1 Mq. 
This range was chosen to help ensure complete sampling 
in all bins since the lowest MSTO mass in our sample cor- 
responds to « 0.75 Mq, and the photometric errors remain 
small (< 0.05 mag) within the magnitude range for each stel- 
lar mass bin in every cluster. We obtained number counts 
for all stellar mass bins in the annulus r c < r < 2r c , where 
r is the distance from the cluster center. This reduced our 
sample size by five clusters since the spatial coverage offered 
by the ACS field of view is incomplete in these cases. 

We obtained completeness corrections for each stellar 
mass bin in the annulus immediately outside the core (r c < 
r < 2r c ). This was done using the results of artificial star 
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Table 2. Completeness-corrected number counts for all five stel- 
lar mass bins in the annulus r c < r < 2r c 



Cluster 


MSI 


MS2 


MS3 


MS4 


MS5 


ID 


(0.65 - 0.75 Mq) 


(0.55 - 0.65 Mq) 


(0.45 - 0.55 M ) 


(0.35 - 0.45 M ) 


(0.25 - 0.35 Mq) 


104 


15113 


16056 


- 


- 


- 


1261 


4652 


4840 


4747 


4647 


4441 


2298 


1022 


865 


828 


673 


595 


4147 


525 


355 


257 


154 


100 


4590 


2830 


3179 


3793 


4477 


5543 


5024 


9773 


9626 


9725 


8114 


5859 


5272 


9763 


10447 


12097 


12886 


14339 


5286 


9394 


8330 


6436 


2765 


941 


5904 


5382 


6726 


9226 


- 


- 


5927 


4304 


4208 


5244 


6043 


- 


5986 


10328 


10936 


12070 


13012 


13635 


6093 


4356 


2658 


- 


- 


- 


6121 


1111 


879 


- 


- 


- 


6171 


1207 


1049 


968 


1064 


- 


6205 


11757 


13012 


16176 




- 


DZlo 




OQQ7 

Zoo ( 


Zo4o 


zDoy 




6254 


4631 


4826 


5375 


6394 


7342 


6304 


1806 


1941 








6341 


4127 


4019 


3771 


2456 




6535 


292 


200 


188 


143 


134 


6584 


3083 


3330 


3624 


3880 


4498 


6637 


3166 


3192 


3710 


3818 


1714 


6779 


3273 


2941 


2983 


3047 


3088 


6838 


592 


634 


654 






6934 


2615 


2416 


2228 


1836 


1302 


6981 


2731 


2774 


2914 


2865 


2847 


7089 


12549 


12388 


10113 


3458 





tests taken from I Anderson e~ al. (2008). Number counts for 
each mass bin were then multiplied by their corresponding 
completeness corrections. We did not include core number 
counts in our analysis since our completeness corrections 
begin to exceed 50% somewhere inside the core for every 
cluster in our sample. This is due to crowding and the high 
central surface brightnesses at the centres of our clusters. 
We have entire ly removed three clusters fro m our original 
sample used in iLeigh. Sills fc Kniggel l|201ll ). namely NGC 
1851, NGC 5139, and NGC 6652. This is because their com- 
pleteness corrections exceeded 50% in every mass bin in the 
annulus immediately outside the core. We also removed ad- 
ditional clusters from our samples for the lowest three mass 
bins whenever their completeness corrections exceeded 50%. 
These clusters typically had the highest MSTO masses. In 
total, this left us with 27, 27, 23, 20, and 15 clusters in each 
of the five mass bins, in order of decreasing stellar mass. The 
completeness-corrected number counts for each stellar mass 
bin have been provided in Table [2] 
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The field of view of the ACS images is about 200" on 
a side, which gives physical scales ranging between 1.5 and 
16 pc (for the closest and furthest clusters in our sample). 
Based on this, we expect foreground contamination by field 
stars to be negligible for most of the clusters in our sam- 
ple given thei r curre nt locations in the Galaxy. For exam- 
ple, [5aCosta| |l98i) considered star count data in a similar 
area and over a comparable range of stellar masses for three 
nearby globular clusters. The author found that the cor- 
rections resulting from field contamination were always less 
than 10% over nearly the entire range of stellar masses we 
are considering. 



2.3 Weighted Lines of Best-Fit 

In order to quantify cluster-to-cluster differences in the 
present-day stellar mass functions of the clusters in our sam- 
ple, we obtained lines of best-fit for (the logarithm of) the 
number of stars belonging to each stellar mass bin versus 
(the logarithm of) the total number of stars spanning all 
five mass bins, which provides a proxy for the total cluster 
mass. This can be written: 



Table 3. Initial Model Parameters 



log N biri 



- Ji log 



Ntot 
10 3 



+ 6 U 



(1) 



where Nbin,i is the number of stars belonging to mass bin i, 
Ntot is the total number of stars spanning all five mass bins, 
and ji and Si are both constants. 

Our motivation for adopting this technique is as follows. 
If the fraction of stars belonging to each mass bin, or fbi n ,i = 
N(,i ni i/Ntot, is constant for all cluster masses, then we would 
expect ~Nbin,i to scale linearly with Ntot. Or, equivalently, ji 
~ 1 in Equation Q] However, if there is any systematic de- 
pendence of fbin,i on the total cluster mass, then we should 
find that Nf, in> 3 does not scale linearly with N to t. In log-log 
space, the slope of the line of best-fit for stellar mass bin i 
should be less than unity (i.e. 7, < 1) if hin.i systematically 
decreases with increasing cluster mass. Conversely, we ex- 
pect 7i > 1 if fbin.i systematically increases with increasing 
cluster mass. This means that, for a sample of clusters with 
a wide range of total masses, we expect 7,: < 1 for the high- 
est mass stars and 7, > 1 for the lowest mass stars. This 
is because clusters lose preferentially low-mass stars due to 
two-body relaxation, and this process operates the fastest 
in lower mass clusters. 

Equation [T] quantifies the number of stars belonging to 
each stellar mass bin as a function of the total cluster mass. 
More generally, it provides a means of quantifying cluster-to- 
cluster differences in the stellar mass function as a function 
of both the stellar mass and the total cluster mass. 

The lines of best-fit have been weighted by adopting 
uncertainties for the number of stars in each mass bin us- 
ing Poisson statistics. Uncertainties for the slopes (i.e. for 
7,: in Equation [TJ were found using a bootstrap methodol- 
ogy in which we generated 1,000 fake data sets by randomly 
sampling (with replacement) number counts from the obser- 
vations. We obtained lines of best fit for each fake data set, 
fit a Gaussian to the subsequent distribution and extracted 
its standard deviation. 



Parameter 


Initial Values 


IMF slope (a) 


1.3, 0.4, 0.0, -0.4 


Number of Stars 


le5, 2e5, 4e5, 6e5, 8e5, le6, 2e6, 4e6 


Concentration (Wo) 


5.0, 5.5, 6.0, 6.5, 7.0 


Virial Radius (in pc) 


3, 4, 5 



2.4 Monte Carlo Models 

We have generated 268 Monte Carlo simulations for globu- 
lar cluster evolution spanning a range of initial total num- 
bers of stars, concentrations, virial radii and IMFs. The 
models realistically take into account both stellar and bi- 
nary evolution, and track both short- and long-range grav- 
itational interactions between both single and binary stars. 
Detailed explanations conc erning the development of these 



model s can be found in Joshi. Rasio fc Portegies Zwarj 



20001), IJoshi. Nave fc Rasiol I J200 111), IFregeau et al.l (|2003h . 



IFregeau fc Rasid (120071 ). and lChatteriee et al.1 i|2010T ). 

For every combination of initial concentration (Wo), 
virial radius (i V i r ) and IMF slope (a in Equation [2}, we 
generated a series of models with different initial total num- 
bers of stars (i.e. total cluster masses). We adopted an IMF 
of the form: 



dN . 
— = /3m 
dm 



(2) 



where a and f3 are constants. This was taken from iKroupal 
(|200ll ). who fit a three-part power-law to this function with 
q = 2.3 for 0.50 < m/M < 1.00, a = 1.3 for 0.08 < m/M Q 
< 0.50, and a = 0.3 for 0.01 < m/M© < 0.08. We varied a 
only in the stellar mass range 0.08 < m/M© < 0.50. 

Each model run was evolved for a period of 12 Gyrs, 
which roughly c oincides with the ages of the clusters in 
our s ample (e.g. |Pe Angeli et al.l l2005: Ma rin-Franch et al.l 
2009). The resulting collection of simulations spanned 
roughly the same range of total masses as our observed sam- 
ple. The initial cluster parameters considered in this paper 
are shown in Table [3] With this suite of simulations, we 
have only scratched the surface in terms of exploring the 
total initial parameter space that could be relevant to the 
GCs in our observed sample. However, our goal in this paper 
is to demonstrate the strength of our technique for quanti- 
fying cluster-to-cluster differences in the observed present- 
day MFs, and to show by example how our method can be 
used to compare a large number of observed MFs to analo- 
gous samples of simulated MFs. We defer a more complete 
exploration of the total possible parameter space of initial 
conditions to a future paper. 

The simulated clusters were placed on circular orbits at 
a distance of 4 kpc from the Galactic centre, and the re- 
sulting tidal effects from the Galaxy were accounted for. We 
note that these effects are reduced by adopting a smaller 
initial virial radius. The effects of tides were typically small 
in all but those models for which both the initial mass and 
concentration were very low, which agrees with th e results 
of previous studies (e.g. IVesperini fe Heggieiri997l ). We as- 
sumed a metallicity of Z = 0.001 for all simulated clusters. 
This roughly agrees with what is typically observed in Galac- 
tic GCs (e.g. lHarrislll996. 2010 update!) , and primarily af- 
fects the rate of stellar mass loss due to winds early on in 
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the cluster lifetime when m assive stars are still present (e.g. 
IChernoff fc Weinberg! Il99ch . Although this does affect the 
rate of dynamical evolution, the effect should be very sim- 
ilar from cluster-to-cluster for the initial conditions consid- 
ered in this paper. This is because the mass loss that occurs 
early on in the cluster lifetime due to stellar evolution has 
the greatest impact on clusters with very low initial con- 
centrations, and can significantly red uce the time for cluster 
dissolu tion in these cases. However, iBaumgardt fc Makind 
l|2003l ) showed that this effect is not severe for the range 
of initial concentrations (Wo = 5-7) considered here. Fi- 
nally, we assumed an initial global binary fraction of 10% 
for all model runs using th e same binary orbital p arameter 
distributions as adopted in lChatteriee et all (|201dh . We will 
return to these assumptions in Section [3] 

We generated simulated CMDs for every model 
run by converting the bolometric luminosity of ev- 
ery star to its corresponding magnitude in the ACS 
F8f4W band. This was done u sing the colour 
sion routine of iPols et al.l 



conver- 

(1998 ), which us e s the spec- 



tral libraries of ILeieune. Cuisinier fc Buserl l|l997l ) and 



iLeieune. Cuisinier fc Buser ( 19981 ). For binary stars, the 
magnitudes of the components were combined in order to 
position them in the CMD as single objects. 

Observations of star clusters are projected onto the 
plane of the sky, whereas the output from our models pro- 
vides only a 3-D distance from the cluster centre for every 
single and binary star. Therefore, it was necessary to convert 
these 3-D distances to corresponding 2-D values. This was 
done by randomly varying the component along the line-of- 
sight to the cluster, and using the 3-D distance to calculate 
a 2-D value. Using these projected 2-D distances from the 
cluster centre, we also generated surface brightness profiles 
for every model run and re-calculated a 2-D core radius (de- 
fined as the distance from the cluster centre at which the 
surface brightness falls to half its central value). These 2-D 
core radii were then used to count the number of objects 
(i.e. single and binary stars) belonging to each stellar mass 
bin located within the annulus immediately outside the core 
(i.e. r c < r < 2r c ). 



o 




log (N r 



< r < 2r ■ 



Figure 1. The logarithm of the number of stars belonging to 
each stellar mass bin (N;,j n ) as a function of the logarithm of the 
total number of stars spanning all five mass bins in the annulus 
immediately outside the core (N r< .< r< 2 r<: ). In descending order 
from top to bottom, the plots correspond to number counts in 
the mass ranges 0.65 - 0.75 M Q (MSI), 0.55 - 0.65 M Q (MS2), 
0.45 - 0.55 M (MS3), 0.35 - 0.45 M Q (MS4), and 0.25 - 0.35 
Mq (MS5). Lines of best fit are shown for each mass bin by solid 
lines. 



Table 4. Lines of Best Fit for log N;, in 4 = (7 ± A7)log 
(N tot /10 3 ) + (S ± AS) 



Stellar Mass 



7 ± A7; 8 ± AS 



MSI (0.65-0.75 M ) 
MS2 (0.55-0.65 M Q ) 
MS3 (0.45-0.55 M©) 
MS4 (0.35-0.45 Mq) 
MS5 (0.25-0.35 Mq) 



0.81 ± 0.09; 2.64 ± 0.12 
0.91 ± 0.09; 2.50 ± 0.10 
0.93 ± 0.03; 2.43 ± 0.04 
0.99 ± 0.05; 2.33 ± 0.08 
1.12 ± 0.11; 2.14 ± 0.17 



2.5 Comparing the Observed and Simulated 
Present-Day Mass Functions 

For every combination of initial concentration, virial radius, 
and IMF slope, the different runs corresponding to different 
initial numbers of stars were grouped together. This gave us 
a sample of MFs spanning a range of total cluster masses for 
every combination of initial conditions. The selection crite- 
ria described in Section ^. 2l was then applied to each model, 
and lines of best-fit were found for each stellar mass. The 
slopes of these lines of best-fit (i.e. ji in Equation [TJ were 
then compared to the corresponding observed slopes for ev- 
ery stellar mass, and both the chi-squared value and the 
probability that the two samples (i.e. the observed 7i's for 
all five mass bins and a given set of theoretically-derived 7i's 
for all mass bins) are drawn from the same distribution were 
found. 



3 RESULTS 

In this section, we present the results of both our observa- 
tional analysis and its comparison to the models. 

3.1 Observational Analysis 

We have plotted the logarithm of the number of stars in each 
stellar mass bin versus the logarithm of the total number of 
stars spanning all five mass bins in Figure [1] The slopes 
and y-intercepts for the weighted lines of best-fit performed 
for each of these relations provided values for 7; and Si in 
Equation [1] These are shown in Table 21 along with their 
corresponding uncertainties (A7 and AS). Each table entry 
has been provided in the form (7 ± A7; 5 ± AS). The values 
for ji have also been plotted in Figure [2] 

As shown in Figure[2] 7* tends to systematically increase 
with decreasing stellar mass. The uncertainty for 7^ is the 



© 2010 RAS, MNRAS OOO.riHTTI 



Quantifying the Universality of the Stellar IMF 7 



T 



1.2 - 
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Figure 2. Slopes for the lines of best-fit, given by 7, in Equa- 
tion [l] plotted as a function of stellar mass. 

highest for the lowest mass bin (MS5 in Table [4]). This is 
because the photometric errors are the highest at these dim 
magnitudes. However, the errors are consistently at most ~ 
10% of the width in magnitude of their corresponding mass 
bin. 

In an attempt to improve upon these statistics, we have 
also calculated reduced chi-squared values with added in- 
trinsic dispersion for the relations for each mass bin. That 
is, for each mass bin we added a constant term to the un- 
certainty for each data point, found the uncertainty that 
yielded a reduced chi-squared of one, and looked at the sub- 
sequent effects on the uncertainties for the line of best-fit. 
Based on this, we appear to be slightly over-estimating the 
uncertainties for the MSI, MS2 and MS3 mass bins using our 
bootstrap approach, and slightly under-estimating them for 
the MS4 and MS5 bins. 

The change in the distribution of stellar masses as a 
function of the total cluster mass can be illustrated using 
pie charts, as shown in Figure [3] Using the values for 7; and 
Si provided in TableSl we have generated pie charts for three 
total numbers of stars (spanning all five mass bins), namely 
N to i = 10 5 , 10 4 , 10 3 (from top to bottom in Figure [3]). As is 
clear, low-mass stars become more and more preferentially 
depleted with decreasing total cluster mass. 

From top to bottom, the pie charts can be interpreted 
as depicting the evolution of the stellar mass function with 
increasing dynamical age. This can be understood as follows. 
The inverse of the half-mass relaxation time can be used as 
a proxy for the rate of two-body relaxation throughout the 
entire cluster. The half-mass relaxation time ranges from 
several million years to the age of th e Universe or longer, 
and is approximated by (|Spitzerlll987l ) : 

t rh = 1.7 x W 5 [r h (pc)f /2 N 1/2 [m/MQr 1/2 years, (3) 
where is the half-mass radius (i.e. the radius enclosing 



half the mass of the cluster), N is the total number of stars 
within Th and m is the average stellar mass. Simulations 
have shown that th changes by a fac tor of at most a few 
over the course of a cluster's lifetime jHenonll 19731 ; iMurravl 
2009). The GCs that comprise the ACS sample show a 
range of masses spanning roughly 3 orders of magnitude 
(10 4 -10 ^ M w ), and have comparably old ages 10-1 2 
Gyrs) (|De Angeli et al.1 120051 ; iMarin-Franch et alj 120091 '). 
Moreover, their half-mass radii typically differ by less 
than a factor of 2 (see Table [TJ. Therefore, Equation [3] 
suggests that the total cluster mass provides a rough proxy 
for the degree of dynamical evolution due to two-body 
relaxation. In other words, the effects of two-body relax- 
ation on the evolution of the stellar mass function should 
be the most pronounced i n the least massive clusters 



in the ACS sample (e.g. Dc Marchi, Paresce & Pulond 


2007: Baume 


ardt, De Marchi & Kroupa 20081; 


Kruiissen & Mieske 


|2009f). Said another wav. dvnami- 



cal age increases with decreasing cluster mass. Therefore, 
Figure [3] shows that our results are consistent with the 
general picture that two-body relaxation is the cause of the 
observed depletion of low-mass stars in low-mass clusters, 
as opposed to some unknown feature of the star formation 
process. 



3.2 Theoretical Analysis 

In this section, we compare the results of our observational 
analysis to 268 Monte Carlo simulations for globular cluster 
evolution spanning a range of initial conditions. 

Figure [4] shows a comparison between all five 7* values 
found from the observed MFs and the corresponding model 
ji values for every combination of initial conditions. As is 
clear, the agreement is excellent for nearly every combina- 
tion of initial IMF, concentration and virial radius. This 
was confirmed by our chi-squared values, and the probabil- 
ity that the observed and model 7;'s are drawn from the 
same distribution exceeded 64% for all comparisons. 

Our results suggest that a Kroupa IMF (i.e. a = 1.3 
in Equation [2} typically gives the best agreement with the 
observations. Every set of models with this IMF yielded a 
probability greater than 93% that the observed and model 
7i's are drawn from the same distribution. Our results also 
appear to be relatively insensitive to the initial concentra- 
tio n and virial radius. This a grees with what was found 
by | Vesperini fc Heggid (|l997h and iBaumgardt fc Makinol 
I2003TI given the limited ranges we have explored for these 
parameters. 



4 SUMMARY & DISCUSSION 

In this paper, we have presented a new technique to quan- 
tify cluster-to-cluster variations in a large sample of observed 
stellar mass functions. Our method quantifies these differ- 
ences as a function of both the stellar mass and the total 
cluster mass, and offers the advantage that it is insensitive 
to the exact functional form of the MF. We have applied 
our technique to completeness-corrected stellar mass func- 
tions in the range 0.25-0.75 Mq for a sample of 27 globular 
clusters taken from the ACS Survey for GCs, and have com- 
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Figure 3. Stellar mass functions depicted in pie chart form. 
The total area of each circle corresponds to the total number of 
stars spanning all five stellar mass bins, and each pie slice shows 
the fraction of this total corresponding to each mass bin. Each 
of these fractions was calculated using the weighted lines of best- 
fit provided in Table [4] From top to bottom, the total number 
of stars used to generate each pie chart was 10 s , 10 4 , and 10 3 . 
Darker pie slices correspond to more massive bins. The sequence 
of pie charts progressing from top to bottom effectively shows the 
evolution of the stellar mass function with increasing dynamical 
age. 
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Figure 4. Comparison of the observed 7; values to the corre- 
sponding model 7^ values for every combination of initial condi- 
tions. The observed 7;'s are shown as solid circles, whereas the 
model 7i's are shown as different symbols. Each inset shows the 
model 7;'s for a different IMF. Starting at the upper right and 
rotating clockwise, the insets correspond to IMF slopes (in Equa- 
tion fSJ of -0.4, 0.0, 0.4, and 1.3. The online version of this plot 
shows the model 7;'s either as open or coloured symbols which 
have been used to indicate different initial concentrations. The 
solid blue circles, solid red squares, solid green triangles, open 
black triangles, and open five-point stars correspond to initial 
concentrations (Wo) of 5.0, 5.5, 6.0, 6.5, and 7.0, respectively. A 
small offset has been implemented for every stellar mass bin to 
indicate different initial virial radii, with the virial radius increas- 
ing from left to right. A horizontal line has also been included on 
every observed data point in order to indicate the width of the 
corresponding stellar mass bin. 



pared the results to a series of Monte Carlo models for GC 
evolution. 

In the subsequent sections, we discuss the implications 
of our results for the formation and evolution of Milky Way 
GCs. 



4.1 The Effects of Two-Body Relaxation 

We have shown that the observed differences in the present- 
day MFs in our sample can be reproduced by assuming (1) 
an universal initial mass function for all clusters, and (2) 
that internal two-body relaxation is the dominant mecha- 
nism contributing to the cluster-to-cluster variations. 

Our results are the most reliable in the mass range 0.45 - 
0.75 Mq. This is due to the larger photometric errors at the 
fainter magnitudes corresponding to lower stellar masses, 
and incompleteness resulting from crowding. Despite the 
high quality of the data used in this study, these issues are 
currently unavoidable given the nature of the observations. 
This will be a key challenge for future studies to resolve, 



however the method we have presented in this paper offers 
a robust means of performing future analyses. We note that 
the resu lts of our observational analysis are consisten t with 
those of lDe Marchi. Paresce fc Portegies Zwartl ((2010) , who 
fit a tapered power-law distribution function with an expo- 
nential truncation to the stellar mass functions of a sample 
of 30 clusters containing both young and old members. We 
h ave also verified tha t our results are consistent with those 
of IPaust et all (|2010h . who performed power-law fits to the 
MFs of 17 GCs taken from the ACS Survey. 

In this paper, we have focused on the local MFs in the 
central cluster regions of the GCs in our sample. However, 
previous studies have shown that the local MF can differ 
consi derably from the global MF (e.g. IVesperini fc Heggiel 
Il997h . In principle, this should not have affected our com- 
parisons between the observed and simulated MFs since we 
have considered the same structural areas of the clusters in 
all cases. 
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4.2 The Effects of Binary Stars 

Binary stars are unresolved in GCs, appearing as single 
objects located above the MS in the cluster CMD. There- 
fore, we have included some objects in our number counts 
that are in fact binaries masquerading as single stars. Pre- 
vious studies have shown that unresolved binaries can con- 
tribute to flattening, or even inverting, the stellar mass func- 
tion in the range of stellar masses considered in this paper 
(e.g. IMarks. Kroupa fc Baumgardtl 120081 ; IMarks fc Kroupal 
2010). Moreover, observational evidence suggests that the 
binary fractions in G Cs are inversely proportional to the to- 
tal cluster mass fe.g.lSollima et al.ll200Sl ; lMilone et al.ll200Sl ; 
iKnigge. Leigh fc Sillsll2009l ) . In particular, the most massive 
MW GCs tend to h ave binary fractions on the order of only 
a few percent (e.g. iRubenstein fc Bailvnl 1 19971 ; ICool et"all 
l2002l ; lDavis et al.ll200Sl ). whereas the least massive MW GCs 
tend to have larger bin ary fractions that c an even exceed 
50% in some cases (e.g. iMilone et al1l201ll ). This suggests 
that unresolved binaries should have had the largest effect 
on the MFs of the lowest mass clusters in our sample. There- 
fore, unresolved binaries could also be contributing to the 
general trend we have found of increasing yt with decreasing 
stellar mass. 

In an effort to quantify the effects of unresolved bina- 
ries on our results, we removed all binaries from our simu- 
lated MFs and re-performed our weighted lines of best-fit. 
This confirmed that unresolved binaries could indeed be con- 
tributing to the general trend we have found of increasing 
7i with decreasing stellar mass. However, this effect was not 
significant in the models, in large part because we assumed 
a constant initial global binary fraction of 10% for all clus- 
ters. In order to reproduce the observed trend in ji, it would 
require that the range in binary fractions between low- and 
high-mass cluste rs is greater than t he observed range by a 
factor > 2 (e.g. IMilone et ail 1201 ll ). Although the binary 
fractions gen erally increase in the core over tim e in our sim- 
ulations fsee iFregeau. Ivanova fc Rasid l|2009h for more de- 
tails) , they are not sufficiently high to have significantly con- 
tributed to the trend of increasing ji with decreasing stellar 
mass. At the end of our simulations, the core binary frac- 
tions are typically in the range 10-30%. 

In order to properly assess these effects, objects that 
are in fact unresolved binaries should be identified and our 
analysis of the observations should be re-performed. This 
could be done using multi-band photometry since, if a given 
binary happens to fall on a single star evolution track in one 
CMD, it is unlikely to fall on the corresponding tracks in 
other CMDs constructed using different wavelength bands. 
Stellar evolution models could then be used to constrain the 
masses of the component stars. We intend to address this 
issue in future work. 

We expect that the influence of binaries on GC evolu- 
tion by, for example acting as h eat-sources via hard e ning 
encounters with single stars (e.g. iHut fc Ba hcall 1983]; iHutl 
1 19831 ; IFregeau. Ivanova fc RasioT 12009 ). should have had a 
negligible impact on our results. This is because most of the 
clusters in ou r sample should still be unde rgoing core con- 
traction (e.g. IGieles. Heggie fc Zhaolfeoill ). It follows that 
their central densities have not yet become sufficiently high 
for encounters involving binaries to occur frequently enough 
that they could have significantly affected the cluster evo- 



lution. Notwithstanding, future studies should incorporate 
models spanning a range of realistic initial binary fractions 
and distributions of orbital parameters in order to properly 
assess all of these effects. 



4.3 The Effects of Tides from the Galaxy 

Tides from the Galaxy effectively reduce the time-scale 
on w hich two-body relaxation operates (e.g. iHeggie fc Hutl 
2003). This primarily serves to make clusters appear more 
dynamically evolved than they otherwise would. The same 
effect is caused by disc shocking which, as with tidal ef- 
fects, should most severely affect clusters with the lowest 
masses and the smallest Galactocentric distances. Therefore, 
the locations of clusters within their host galaxies has been 
shown to play an important role in det ermining the degree 
of fla ttening of their stellar MFs (e.g. IVesperini fc Heggiel 
119971 ). In an effort to quantify the effects caused by tides 
on our observational results, we performed several cuts 
in perigalacticon distance and re-performed our weighted 
lines of best-fit. Esti mates for the perigalacticon di s tances 
were obtained from Dinescu. Girard fc van Altenal (|l999h 
and ICasetti-Dinescu et alj l|2007l) for every cluster in our 
sample. Despite removing clusters with small perigalacticon 
distances for which it is typicall y argued that tidal e ffects 
should be the most severe (e.g. IHeggie fc Hutll2003l ). our 
lines of best-fit remained the same. We caution that these 
effects have not been fully accounted for in the simulations of 
GC evolution performed in this paper, for which we adopted 
circular orbits and a Galactocentric distance of 4 kpc in all 
cases. We intend to adopt a more realistic distribution of 
Galactic orbits for our model clusters in a future paper in 
order to properly assess the effects caused by tides in de- 
termining the present-day MFs of the GCs in our observed 
sample. 



4.4 The Effects of Primordial Gas Expulsion 

As discussed in Marks. Kroupa fc Baumgardtl l|2008l ) and 
IMarks fc Kroupa ( 2010l ). the expulsion of primordial gas 
early on in the cluster lifetime can have a dramatic effect 
on the stellar MFs of clusters. In particular, these authors 
showed that clusters that began their lives with smaller con- 
centrations are more likely to have lost a larger fraction of 
their low-mass stars as a result of this effect. Primordial 
gas expulsion could therefore contribute to improving the 
agreement we have found between the observed and sim- 
ulated MFs. This would be accomplished if our 7, values 
for the models were simultaneously increased for low stellar 
masses and decreased for high stellar masses, as is evident 
from Figure [4] This would occur if primordial gas expul- 
sion had a more pronounced effect on the MFs of the lowest 
mass clusters in our sample. This is not unreasonable since 
the depth of the cluster potential increases with increasing 
cluster mass. 



4.5 Future Work 

Given the very old ages and therefore low metallicities 
([Fe/H] « -2.28 - (-0.37)) of the clusters in our sample, our 
technique could potentially be used to better constrain the 



© 2010 RAS, MNRAS OOO.rTlfTTl 



10 Leigh et al. 



IMFs of old massive star clusters and, more generally, star 
formation in the very early Universe. This will be done in 
a forthcoming paper by considering a larger range of IMFs 
and initial cluster conditions than we have considered here. 

Finally, we wish to point out that the method we have 
presented can be generalized to compare any large sample of 
distribution functions. We intend to illustrate this in a future 
paper by using our technique to quantify cluster-to-cluster 
differences in the orbital distributions (period, eccentricity, 
and mass-ratio) of the binary populations in GCs as a func- 
tion of the total cluster mass. 
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